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Abstract 

Charmonium production at heavy-ion colliders is considered within the comovers interaction 
model. The formalism is extended by including possible secondary J / 1\) production through recom- 
bination and an estimate of recombination effects is made with no free parameters involved. The 
comovers interaction model also includes a comprehensive treatment of initial-state nuclear effects, 
which are discussed in the context of such high energies. With these tools, the model properly de- 
scribes the centrality and the rapidity dependence of experimental data at RHIC energy, y/s = 200 
GeV, for both Au+Au and Cu+Cu collisions. Predictions for LHC, ^fs = 5.5 TeV, are presented 
and the assumptions and extrapolations involved are discussed. 
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1 Introduction 



Disentangling effects related to the production of charmonium in hadronic collisions has been a major 
task for both experimentalists and theoreticians for the last two decades. The discovery of J/ip 
suppression, with respect to heavy lepton pair production, in proton-nucleus (pA) collisions has been 
interpreted as a result of the multiple scattering of a cc pair escaping the nuclear environment - 
the so-called nuclear absorption. Moreover, in high energy nucleus-nucleus (^4+^4) collisions one 
hopes to achieve such large temperatures and densities that a new state of deconfined QCD matter, 
the quark-gluon plasma (QGP), is produced. It was suggested that, in the presence of a QGP, 
the charmonium yield would be further suppressed due to color Debye screening [I]. Indeed, such 
an anomalous, compared to absorption, suppression was observed in Pb+Pb collisions at top SPS 
energy [2]. Alternatively, the SPS experimental results can also be described in terms of final state 
interactions of the cc pairs with the dense medium created in the collision, the so-called comovers 
interaction model (CIM). This model does not assume thermal equilibrium and, thus, does not use 
thermodynamical concepts. Within this model the SPS experimental data can be reproduced with an 
effective dissociation cross section a co = 0.65 mb [3]. 

The theoretical extrapolations to collider energies are guided mainly by two trends. On the one 
hand, models assuming a deconfined phase during the collision pointed to the growing importance of 
secondary J/tp production due to recombination of cc pairs in the plasma. The total amount of cc 
pairs is assumed to be created in hard interactions during the early stages of the collision. Then, either 
using kinetic theory and solving rate equations for the subsequent dissociation and recombination of 
charmonium [U [5], or assuming statistical coalescence at freeze-out [6l El E]> one obtains the final 
J/ip yield. These models predict a disappearance of the charmonium suppression with rising collision 
energy as the lifetime of the plasma phase is expected to grow accordingly. See, however, [9] where 
thermal equilibrium is not assumed and charmonium dissociation is quite large at RHIC. On the other 
hand, the CIM with only dissociation of J/tp's predicts [TO] a stronger suppression at RHIC than at 
SPS due to a larger density of produced soft particles in the collision. It also predicts a stronger 
suppression at y = (where the comovers density is maximal) than at forward rapidities. 

In this context, measurements of J/ip production in Au+Au collisions at RHIC, ^fs = 200 GeV, 
gave interesting although surprising results - the suppression at mid-rapidity was on the same level as 
at SPS [lit I12j. This was also the case for Cu+Cu collisions at the same energy |13j . Furthermore, 
the suppression at forward rapidity in Au+Au collisions was stronger than at mid-rapidity. The latter 
feature was not seen for the much smaller collision system created in Cu+Cu collisions. 

The CIM is based on the well known gain and loss differential equations in transport theory. The 
introduction of a recombination term is actually required for detailed balance. So far it had not been 
introduced in the model just because its effect was assumed to be small. The aforementioned RHIC 
results prompt us to a careful evaluation of the effect of the recombination (gain) term. If its effect is 
not negligible it will increase the final J/ip yield. Moreover, this increase will be larger at y = than 
at forward rapidities due to the narrow rapidity distribution of charm, and thus complies with RHIC 
data. 

In the present work we will extend the CIM by allowing recombination of cc pairs into secondary 
J/^'s. We will estimate this effect using the density of charm in proton-proton (pp) collisions at 
the same energy and at various rapidities. Therefore, the model does not involve any additional 
parameters. In Section [2] we present the details of the model; both nuclear effects related to the initial 
state at high energy and final state interactions are described. We calculate recombination effects 
at mid- and forward rapidities in Au+Au and Cu+Cu collisions at RHIC energy y/s = 200 GeV in 
Section [3J The suppression, found without any tuning of the parameters in the model, is in good 
agreement with the data. In Section 0] we make predictions for Pb+Pb collisions at LHC, y^s = 5.5 
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TeV, and discuss uncertainties related to this extrapolation. Finally, conclusions and final remarks 
are given in Section 



2 Description of the model 

In this section, we give a short and updated description of the CIM which has been used to make 
predictions for J/ tp production in A +A collisions at RHIC [10] . This model contains a comprehensive 
treatment of initial-state nuclear effects, such as nuclear shadowing and nuclear absorption, and final 
state interactions with the co-moving matter. Here we will extend the model as follows 

• update nuclear shadowing for hard production of charmonium calculated from parameterization 
of diffractive gluon distribution function [T3] 

• extend nuclear absorption and its energy dependence to the whole kinematically allowed region 

• include the possible recombination of cc into secondary J/tp's, i.e. the gain term in the differential 
rate equation of the model. 

The latter effect is the main novel feature that will enable us to calculate the suppression pattern 
at RHIC, and make further predictions for the upcoming heavy-ion runs at LHC. Nuclear effects 
in nucleus-nucleus collisions are usually expressed through the so-called nuclear modification factor, 
R J Jg(b), defined as the ratio of the J/ip yield in A+A and pp scaled by the number of binary nucleon- 
nucleon collisions, n(b). We have then 



n(b) dN^/dy 
J d 2 s a AB {b) n(b, s) Sf^jb, s) S abs (b, s) S co (b, s) 
J d 2 sa AB (b)n(b,s) 

2, 



(1) 



where a A B(b) = 1 — exp[— a pp AB Tab(6)], Tab{V) = f d sT A {s)TB{b — s) is the nuclear overlapping 
function and T A (b) is obtained from Woods-Saxon nuclear densities [15] . In eq. ([T|), 

n(b,s) = a pp ABT A (s)T B (b- s)/a AB (b) , (2) 

where upon integration over d 2 s we obtain the number of binary nucleon-nucleon collisions at impact 
parameter b, n(b). 

The three additional factors in the numerator of eq. (pQ), S sh , S abs and S co , denote the effects of 
shadowing, nuclear absorption (both initial-state effects) and interaction with the co-moving matter 
(final-state effect), respectively. They will be defined below. 



2.1 Initial-state nuclear effects 

The effects related to particle production in hadronic collisions off nuclei are often called initial-state 
effects and have been extensively treated in the literature in the case of J ftp production at different 
energies [16 1 117 } [T8]. Since the origin and relevance of these effects is still under debate, we will present 
a short and updated discussion of them at both low and ultra-relativistic energies. 

At low energies the primordial spectrum of particles created in scattering off a nucleus is mainly 
altered by (i) interactions with the nuclear matter they traverse on the way out to the detector and 
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(ii) energy-momentum conservation. The first effect is called nuclear absorption and is usually pa- 
rameterized within a probabilistic Glauber model. The latter is inferred from the rapidity dependence 
of the spectra. For A+A collisions, these effects can be combined into the following suppression factor 

q abs = [1 - exp(-C(x + )a QQ AT A (b))] [1 - exp(-^(x-)a QQ BT B (b - s))] 

t{x + )^)a 2 Q ABT A {s)T B {b-s) ' U 

where S,(x±) = (1 — e) + ex± determines both absorption and energy-momentum conservation^ In 
[T6] it has been found that 7 = 2, e = 0.75 and o q q = 20 mb give a good description of data. This 
corresponds to a a b s = 5 mb at mid-rapidity which is in agreement with other studies. 

Several novel features are expected to appear at high energies due to the change in the space-time 
picture of particle production. First of all, the form of eq. ([3]) will change due to coherence effects 
|19j . Motivated by results from deuteron-gold collisions at y/s = 200 GeV |20j . it has recently been 
realized that this change should be accompanied by a transformation of £(x±)ctqq into ex\a Q Q. In 
this way a a bs = at high energies [21] [22]. Nevertheless, S abs 7^ 1 for x± not too small, due to 
energy-momentum conservation. This leads to a model for the high-energy transition in which the 
factors in the numerator of eq. ([3]) are changed according to 

I _ e -^±)a Q QAT A {b) 

— > AT A (b) e W [-ea QQ xlAT A (b)/2] . (4) 

t,{X±)<J Q Q 

In the central rapidity region at RHIC, x + ~ 0.025 — 0.05 and absorption effects can be discarded. 
At \y\ > 0, S abs / 1 due to energy-momentum conservation. This effect is usually neglected in the 
discussion of A+A collisions. 

Secondly, coherence effects will lead to nuclear shadowing for both soft and hard processes at 
RHIC, and therefore also for the production of heavy flavor. Shadowing can be calculated within 
the Glauber-Gribov theory [23J, and we will utilize the generalized Schwimmer model of multiple 
scattering [21]. In this case the second suppression factor in eq. ([T]) is given by 

= 1 + AF(y A )T A (s) 1 + BF(y B )T B (b-s) ' (5) 

where the function F{y) encodes the dynamics of shadowing and will be discussed in detail below. 
For a general discussion of nuclear shadowing, see \2b\ 126] . 

Parameterizations of diffractive structure function from 7*iV scattering as measured at HERA 
have been utilized in [27] to find shadowing for sea quarks. Most recently gluon shadowing has been 
calculated in [28] using recent data on diffractive gluon density function [14] . Then shadowing is 
governed by 

F(x,Q 2 ) = 4vr / F dxjpBixjp) t2 \ [ ^%f Fl(t min ) , (6) 

where F(x,Q 2 ) is related to F{y) in eq. (|5|) through kinematical relations. For gluon fusion, x = 
rriT exp(±y)/- v /i. We put Xp ax = 0.1, where shadowing is expected to disappear. In eq. ([6]), ^(x, Q 2 ) 
is the structure function for a nucleon, F^(x]p, Q 2 ,[5) is the t-integrated diffractive structure function 
of the nucleon, B(xjp) is the slope parameter of the distribution, and F A {t m i n ) is the nuclear form factor 
where t m i n — m 2 N xjp. Equation ([5]) determines shadowing for quarks (antiquarks) in nuclei [27]. For 
gluons the same expressions have been used with the substitutions: F^ (xjp, Q 2 , 0) — > Fjp (xjp, Q 2 , 0) , 



recall that x± = (y/x% - AM 2 / s ± x F )/2. 
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i*2 (x, Q 2 ) — > xg (x, Q 2 ), Fp and g represent gluon distributions in the Pomeron, measured in diffrac- 
tive deep inelastic scattering, and in the proton, respectively. We use this information to calculate 
shadowing for J ftp assuming that it originates from a color octet cc pair produced in a hard gluon 
fusion process with Q = tut = 4 GeV. The same value of Q is used for open charm production. 

Since experiments at HERA mostly deal with hard diffraction and their parameterizations are 
quite uncertain for Q 2 < 4 GeV 2 , the density of comovers, mostly low-p± particles, will be calculated 
in the spirit of the model presented in |25t 129] . where shadowing corrections are given without free 
parameters in terms of the triple-Pomeron coupling determined from diffraction data. Then, in eq. (|5]) 

F(y A ) = C sh [e W [AY± x ) - exp(AFA n )] , (7) 

where Y^ in = ln(i? J 4mjv/v / 3) ) A = 0.13 and C s h = 0.31 fm 2 . The value of Y max depends on the 
rapidity of the produced particle h, Y^ ax = ln(s/mj.)/2 ± y with the +(— ) sign if h is produced 
in the hemisphere of nucleus B(A). mx is the transverse mass of the produced particle. We use 
rriT = 0.4 GeV at RHIC and 0.5 GeV at LHC. This model has been used to correctly predict the 
centrality dependence of soft particle production at RHIC [29j . 

The energy dependence of these different effects can be summarized as follows. Particle production 
at SPS is dominated by low-energy effects, i.e. nuclear absorption given by eq. (J3j) and small nuclear 
shadowing, while RHIC already belongs to the high-energy regime. Nuclear shadowing is non-negligible 
at mid-rapidity, and the combined effect of shadowing and energy-momentum conservation, in the 
spirit of eq. should be accounted for at forward rapidities. At LHC, shadowing will be very strong 
even at y = while energy-momentum conservation is a small effect and can be neglected in most of 
the kinematics. We will now proceed with the discussion of final state effects. 

2.2 Dissociation by comovers interaction and recombination 

The CIM was developed to explain both the suppression of charmonium yields [3l fT0 | l30~ t l3T j 1321 1331 IM] 
and the strangeness enhancement |35U36| in nucleus-nucleus collisions at the SPS. Neglecting possible 
recombination effects, the rate equation governing the density of charmonium in the final state, Nj/^p, 
can be written in a simple form assuming a pure longitudinal expansion of the system and boost 
invariance. For an A+A collision the density of J/ij) at a given transverse coordinate, s, impact 
parameter b, and rapidity is given by 

T ^JJ± = -a co N c °(b,s,y)N J/4) (b,s,y) , (8) 

where a co is the cross section of charmonium dissociation due to interactions with the co-moving 
medium, with density N co . It is found from fits to low-energy experimental data to be a co = 0.65 mb 
[3]. To incorporate the effects of recombination, we have to include an additional gain term propor- 
tional to the (squared) density of open charm produced in the collision. Then eq. (jSJ) is generalized 
to 



dr 



(b,s,y) = -a co N co (b,s,y)Nj /i; {b,s,y) - N c {b, s,y)N 5 {b, s,y) , (9) 



where we have assumed that the effective recombination cross section is equal to the dissociation cross 
section. 0. This extension of the model therefore does not involve additional parameters. 
Equation ([9j) cannot be solved analytically. We approximate its solution as 



S co (b,s,y) = exp < -a c 



N J/f (b,s,y) 



hi 



N co (b,s,y) 
N pp (0) 



(10) 



These two cross-sections have to be similar but not necessarily equal. We have taken the simplest possibility. 
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resembling the exact solution of eq. ([8]), since the first term in the exponent of eq. (I10|) is exactly the 
survival probability of a J ftp interacting with comovers [10] . The density of open and hidden charm in 
A+A collisions, N c , N 5 and Nj/^, respectively, can be computed from their densities in pp collisions as 
N AA (b,s) = n(b, s)S^Q(b, s)N pp , with similar expression for N AA and NjA. Here n(b,s) is given by 

eq. (j2J) and Sjjq is the shadowing factor for heavy quark production, given by eq. 0. Then eq. (fTUj) 
becomes 



S c °(b,s,y) 



exp 



-Or. 



N co (b,s,y)-C(y)n(b,s)S s H h Q (b, 



where 



C{y) 



(dN^/dyf (da%/dy 



In 



N co (b,s,y) 
N pp (0) 



(11) 



dNp^ I dy o-pp dcr//^ / dy 



(12) 



The quantities in the rightmost term in eq. (|12p are all related to pp collisions at the corresponding 
energy and can be taken from experiment or a model (for extrapolation of the experimental results). 
The cc pairs are mostly in charmed mesons, such as D and D*. With a co fixed from experiments 
at low energy, where recombination effects are negligible, the model, formulated above, should be 
self-consistent at high energies^. We expect the effect of recombination to be stronger at mid- than 
at forward rapidities. At y ^ the recombination term is smaller (relative to the first one) since the 
rapidity distribution of D, D* is narrower than the one of comovers. This will produce a decrease of 
R"aa w ith increasing y which may over-compensate the increase due to a smaller density of comovers 
at y ^ 0. 

The density of comovers is calculated using the dual parton model [37] together with the proper 
shadowing correction. Then 



N co (b,s,y) = N c J s (b,s,y)S$(b, S ,y) , 



(13) 



where denotes the shadowing for light particles and is given by eqs. ([5]) and ([7J). The non-shadowed 
(NS) multiplicity of comovers is then 



Nf?s(b,8,y) 



3 dAffi s (M,y) 

2 dy 



where 



n A (b,s) 



AT A {s) 



{Ci(6,y)n A (M) + C 2 (b,y)n(b,s)} , 



[1 - exp (-a pp BT B (b - s))] 



(14) 



(15) 



and n(b,s) is given by eq. ([2]). The numerical values of the coefficients C\ and C2 can be found in 
[TO] . C\ drops with energy while C2 increases, and at LHC Ci ~ while C2 ~ 6 at mid-rapidity. Also 
in eq. (IQJ) N pp (0) = l{dN ch /dy) p JL /irR 2 p , which we estimatf] to be 2.24 fin" 2 at y/s = 200 GeV and 
4.34 fm~ 2 at v/i = 5.5 TeV. Finally, (Tpp, taken as its non-diffr active value, is 34 mb and 59 mb at 
RHIC and LHC, respectively. 



Note, however, that a co could change when the energy increases. We do not expect this effect to be important and, 
since we are unable to evaluate the magnitude of this eventual change, we shall use the same value a co — 0.65 mb at all 
energies. 

4 pp values at LHC of dN/dy and a are based on DPMJET-III calculations. We thank J. Ranft for providing these 
results. 



6 



3 Charmonium production at RHIC 



At RHIC, the dissociation term alone gives a too strong suppression compared to experimental data 
|10j . We therefore proceed to estimate the effect of recombination. The density of open charm at mid- 
rapidity in pp collisions at y/s = 200 GeV has been reported in [38] and the most recent measurement 
of the J/ip density in [HU]. We present results of experimental measurements in Table [TJ The semi- 
leptonic branching ratio for the J/ip is Bu = 0.059. Then, from eq. (|12f) it follows that C = 0.59 at 
mid-rapidities. One has to keep in mind, however, that the measured yield of open charm at RHIC 
is almost twice as large as predictions from pQCD [40J. In the left picture of Fig. [1] we present the 
results of our model compared to experimental data at mid-rapidity. The different contributions to 
J/ip suppression are shown (see figure caption). Note that at mid-rapidities the initial-state effect is 
just the shadowing. As discussed above nuclear absorption due to energy-momentum conservation is 
present at forward rapidities but is negligibly small at mid-rapidities. 

Measurements of open charm at forward rapidities [4T] have too large systematic errors at the 
moment. Therefore, to estimate the density of charm at forward rapidity we use results of PYTHIA 
|42j with parameters and settings as described in [41] . The resulting rapidity distributions can be 
described by slightly broadened Gaussian functions. Note, that with this PYTHIA value, given in 
Table 1, the ratio of open charm production at mid- and forward rapidities is very similar to the 
one measured for J/ip production. With these values of J ftp and open charm rapidity distributions 
we get C = 0.32 at forward rapidities. Using this value of C and a ratio of 1.2 [43] of the comover 
densities N co between y = and y forward we obtain the curve in the r.h.s. of Fig. 1. We note 
that no free parameters were tuned to obtain these results. Note that, contrary to the results in [10] 
with no recombination, the J/ip suppression at forward rapidity is somewhat larger that the one at 
mid-rapidities, in agreement with experimental data. This is due both to the recombination term and 
to the initial-state effects. The latter are stronger for forward rapidities. They include the effects of 
energy-momentum conservation, which were not considered in |10j . 

For consistency, we have also made calculations for the J/ip suppression in Cu+Cu collisions 
at RHIC which was reported in [13| . using the same parameters as above for Au+Au collisions. The 
results are shown in Fig. [21 and are in good agreement with the experimental data, except for peripheral 
collisions where the error bars are quite large. 

Concluding, our procedure gives a reasonable description of data both at mid- and forward rapidity 
for different collision systems at RHIC, in contrast with the recent claim |44] that the CIM fails to 
reproduce forward rapidity data. We have shown that the situation is much improved by inclusion of 
a recombination term. 

4 Predictions for LHC 

Based on our previous discussion, it is obvious that recombination effects will be of crucial importance 
in Pb+Pb collisions at LHC {y/s = 5.5 TeV). Assuming that the energy dependence of open charm 
and J/ip in pp collisions is the same (between RHIC and LHC energies), the energy dependence of 
the parameter C will be that of crp^,/a P p. The total and differential cross section for charm can be 
calculated using perturbative techniques [401 145] . The calculations for low energies are in agreement 
with data, yet predictions for RHIC and Tevatron energies are lower than the data. Therefore, the 
extrapolation to LHC is quite uncertain. If we parameterize the energy dependence of open charm 
production as a cc oc s a , with a = 0.3 and use the values of non-diffractive a pp given at the end of 
Section 2 we obtain C = 2.5 at LHC - a value about four times larger than the corresponding one at 
RHIC. In view of that we consider that realistic values of C at LHC are of the range 2 to 3. In Fig. [3] 
we have calculated the J/ip suppression at LHC for several values of C, including the case of absence 
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of recombination effects (C = 0). Although the density of charm grows substantially from RHIC to 
LHC, the combined effect of initial-state shadowing and comovers dissociation appears to overcome 
the effect of parton recombination. This is in sharp contrast with the findings of [6], where a strong 
enhancement of the J ftp yield with increasing centrality was predicted. 

It is clear that a better theoretical control is needed on the various factors that are included in C. 
As discussed in Section 2 the comovers cross section a co can vary a little with energy. But the most 
important theoretical input at the moment is the energy dependence of the total charm cross section 
in pp collisions. 

5 Conclusion 

In this work we have incorporated the effects of recombination of cc pairs into J/ip in the comovers 
interaction model. These effects are negligible at low energies (SPS) due to the low density of open 
charm. This model does not assume thermal equilibrium of the matter produced in the collision and 
includes a comprehensive treatment of initial-state effects, such as shadowing, nuclear absorption and 
energy-momentum conservation. 

We estimate the magnitude of the recombination term from J /ip and open charm yields in pp 
collisions at RHIC. Without any adjustable parameters, the centrality and rapidity dependence of 
experimental data is reproduced both for Au+Au and Cu+Cu collisions. 

Finally, we make predictions for future measurements of J/ip in Pb+Pb collisions at LHC. In our 
approach, the magnitude of the recombination effect is controlled by the total charm cross section in 
pp collisions, and therefore our predictions are strongly dependent on input from theoretical models 
at these energies. For a reasonable choice of parameters, we predict that the suppression observed at 
RHIC and lower energies will still dominate over the recombination effects. This is due to the large 
density of comovers and to the strong initial-state suppression at these ultra-relativistic energies. 
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Results for J /if) suppression in Au+Au at RHIC (i/i = 200 GeV) at mid- (left figure), 
and at forward rapidities (right figure). Data are from The solid curves are the 
final results. The dashed-dotted ones are the results without recombination (C = 0). 
The dashed line is the total initial-state effect. The dotted line in the right figure is the 
result of shadowing. In the left figure the last two lines coincide (see main text). ... Il2 
Results for J /if) suppression in Cu+Cu at RHIC (y^s = 200 GeV), at mid- (left) and at 

forward rapidities (right). For details, see caption of Fig. [H Data are from [13] Il3 

Results for J/ip suppression in Pb+Pb at LHC (i/i = 5.5 TeV) at mid-rapidities for 
different values of the parameter C. The upper line is the suppression due to initial-state 
effects (shadowing) [14 
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Figure 1: Results for J ftp suppression in Au+Au at RHIC (yfs = 200 GeV) at mid- (left figure), 
and at forward rapidities (right figure). Data are from [llj. The solid curves are the final results. 
The dashed-dotted ones are the results without recombination (C = 0). The dashed line is the total 
initial-state effect. The dotted line in the right figure is the result of shadowing. In the left figure the 
last two lines coincide (see main text). 
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Figure 2: Results for J/tfi suppression in Cu+Cu at RHIC (y/s = 200 GeV), at mid- (left) and at 
forward rapidities (right). For details, see caption of Fig. [TJ Data are from |13j . 
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Figure 3: Results for J/ip suppression in Pb+Pb at LHC {-\fs = 5.5 TeV) at mid-rapidities for different 
values of the parameter C. The upper line is the suppression due to initial-state effects (shadowing). 
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Table 1: Open charm and J/tfi production cross sections in pp collisions at = 200 GeV. Data are 
taken from [5511551 HI]. 









EXP 


(dafp/dy) 


PYTHTA 




[B a da^/dy) 


EXP 


mid-rap 


123 ± 12 ± 45 jih 




44.3 ± 1.4 ± 5.1 nb 


forward 




70.9 ± 14 


27.61 ± 0.37 ± 0.83 nb 
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